/*
*  Rozenas's Ambiguity Model   
* -----------------------------------------------------
*
*  Author: Lukas Stoetzer  
*  Date: 07 May 2018
*
*/

data {

  // For Observed Data
  int<lower=1> J; // # parties
  int<lower=1> N; // # responses
  int<lower=1> I; // # experts
  int<lower=0> T; //  # categories
  
  int<lower=1,upper=J> jj[N]; // # party for n
  int<lower=1,upper=I> ii[N]; // # expert for n
  
  ordered[T-1] c;             // # threshold paramters; fixed 
  
  int<lower=0> y[N];          // # response on left-right for n; y in {1 ... T}

  
  // For Missig Data
  int<lower=1> N_miss; // # responses
  
  int<lower=1,upper=J> jj_miss[N_miss]; // # party missing for n
  int<lower=1,upper=I> ii_miss[N_miss]; // # expert missing for n
  
  int<lower=0, upper=1> y_miss[N_miss];   // # missing response on left-right for n; y_miss in {0,1}
  
}


parameters {

  vector<lower=0>[I] phi;         // # expert discrimination free
  vector[I] tau;           // # expert difficultly free

  vector[J] mu;                 // # party platform
  vector<lower=0>[J] nu;        // # party ambiguity
  
  vector[2] alpha;

}

// transformed parameters {
// 
//   vector[I] phi;             // # expert discrimination constrained
//   vector[I] tau;             // # expert difficultly constrained
//   tau = append_row(tau_free, rep_vector(-1*sum(tau_free), 1)) ;
//   phi = phi_free * pow(exp(sum(log(phi_free))), (-inv(I)));
// 
// }

model {
  
    // Auxiallary vectors and reals
    vector[T] prob;               // probabailites
    real zs;                      // linear predictor
    real si;                      // denominator of linear predicator
    vector[N_miss] eta;

    // Prior Specifications
    mu ~ normal(0, 1);            // Prior on platform with fixed hyper paramters for identification.
    nu ~ normal(0, 1);            // Prior on ambiguity trunacted at 0.
    
    tau ~ normal(0,1);            // prior on difficultly
    phi  ~ lognormal(0,1);        // prior on discrimination with expecation of 1
    
    alpha ~ normal(0,1);         // prior on irt model for missings

    // Model for positions; loop over all observations
    for (n in 1:N){

      // Model for latent trait
      // zs = mu[jj[n]] ;
      zs = phi[ii[n]]* mu[jj[n]] + tau[ii[n]];
      si = phi[ii[n]]* nu[jj[n]];
      
      // Odered probabilities
      prob[1] = Phi_approx((c[1] - zs)/si);
      for (t in 2:(T-1)){
        prob[t] =  Phi_approx((c[t] - zs)/si) - Phi_approx((c[t-1] - zs)/si);
      }
      prob[T] =  1 - Phi_approx((c[T-1] - zs)/si);
   
      // Observed outcome
      y[n] ~ categorical(prob);

    }
    
    // Model for missigness; loop over all observations
    for (n in 1:N_miss){
      // probit model
      eta[n] = Phi_approx(alpha[1] + alpha[2] * nu[jj_miss[n]] * phi[ii_miss[n]]);
    }
    y_miss ~ bernoulli(eta);


} 
